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Abstract 

We present next-to-next-to-leading order (NNLO) parton distribution functions (PDFs) from 
the CTEQ-TEA group. The CTIONNLO PDF fit is based on essentially the same global data 
sets used in the CTIO and CTIOW NLO PDF analyses. After exploring the goodness of the fits 
to the HERA combined data and the Tevatron jet data, we present various predictions at NNLO 
accuracy for both existing and forthcoming precision measurements from the CERN Large Hadron 
Collider. The range of variations in the gluon distribution introduced by correlated systematic 
effects in inclusive jet production is also examined. 
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I. INTRODUCTION 



A global analysis of perturbative QCD makes use of experimental data from many short- 
distance scattering processes to construct, within some approximation, universal parton 
distribution functions (PDFs) for the proton. These PDFs can be used to calculate hadronic 
cross sections. The CTEQ global analysis at the next-to-leading order (NLO) in the strong 
coupling constant has been developed over decades. Examples of general-purpose PDFs 
in this series include CTEQ6, published in 2002 Q, followed by CTEQ6.1 in 2003 11. 
Starting from CTEQ6.5 published in 2006 [3|, and in ensuing PDF sets such as CTEQ6.6 [J] 
and CT09 jll, the effects of finite quark masses on the CTEQ global analysis have been 
implemented in the S-ACOT-x factorization scheme j6|-0] at NLO accuracy. The most 
recent CTEQ NLO PDFs, named CTIO and CTIOW, were pubhshed in 2010 [lo| and are 
currently in wide use in phenomenological predictions for the Tevatron, LHC, and other 
experiments. 

The data from the CERN Large Hadron Collider (LHC) cover a wide kinematic range with 
high expected precision; a similar level of precision is needed for the theoretical predictions. 
Thus, there is a need for using cross sections calculated up to next-to-next-to-leading order 
(NNLO) in the strong coupling constant with parton distribution function sets that are 
also determined at NNLO. Complete calculations for this order in are available for the 
running coupling aJQ), PDF evolution in Q [ll|, I2I , m atrix elements in deep- inelastic 
scattering (DIS) [l3-19|, and vector boson production jiol. 2ll |. NNLO matrix elements are 
unknown for several other processes in the global analysis. Most notably, the key input 
theoretical cross section for inclusive jet production in pp ' collisions is still evaluated at 
NLO, although some of the NNLO radiative contributions have been already computed 



22|, [231.^ Various NNLO PDF sets have been published in the literature (24l29| that make 



use of NNLO matrix elements when available. 

In this paper, we present a new generation of CTEQ parton distributions, at NNLO, 
named CTIO NNLO.^ They are obtained from a global analysis of QCD data in which three 
changes have been made to include the calculations at NNLO. First, the parton distribution 
functions f{x,Q) are evolved according to the 3-loop Dokshitzer-Gribov-Lipatov-Altarelli- 
Parisi (DGLAP) equations. Second, the strong coupling as(/iij) evolves in the renormaliza- 
tion scale fiR according to the 3-loop beta function. Third, the hard matrix elements for 
DIS and vector boson production are calculated up to two QCD loops. A detailed imple- 
mentation of neutral-current DIS cross sections in the S-ACOT-x scheme, which performs 
a consistent treatment of non-zero masses of heavy partons, up to this order of accuracy is 



documented in Ref. (3l| and recapped in Sec. [TTl 



The experimental data sets included in the CTIONNLO fit are essentially the same as in 
the CTIO and CTIOW NLO fits, with the exceptions of Tevatron Run-1 inclusive jet data 
and a subset of the Tevatron Run-2 lepton charged asymmetry data from W boson decays. 



^ Some other cross sections, such as for Higgs boson production and diboson production, have been calcu- 
lated to NNLO, but to date have not been implemented in any global analysis. 
2 Parametrizations of the CTIONNLO PDFs were publicly released in 2012 [30|. 



3 



The changes in the selection of the experimental data sets are summarized in Sec. IIV[ No 
LHC data are included in the CTIO NNLO analysis, which can therefore be used to make 
predictions based exclusively on the pre-LHC data. [The impact of the LHC data will be 
investigated in the post-CTlO fits.] 

As stated previously, the inclusive jet cross sections in the CTIONNLO analysis are still 
evaluated with NLO matrix elements, hence require special scrutiny. The experimental data 
on the inclusive jet cross sections are statistically very precise but have significant systematic 
uncertainties. We examine variations in the gluon and other PDFs caused by the QCD scale 
dependence and by various treatments of systematic uncertainties in inclusive jet production. 
The impact of jet-related uncertainties is compared to those from other sources. 

For consistency with the CTIO NLO PDFs, CTIO NNLO assumes the same value of the 
QCD coupling strength as{Mz) = 0.118 and pole masses for heavy quarks of rric = 1.3 GeV 
and '^]dL = 4.75 GeV. The CTIONNLO eigenvector PDF sets are available on the CTEQ web 
site jlll, and in the LHAPDF standard format [33|. Together with the PDF eigenvector 
sets for the central value of as{Mz) = 0.118, an additional PDF series in which as{Mz) is 
varied in the range 0.112-0.127 are provided. These PDFs are sufficient for computing the 
correlated PDF-l-a... uncertainty by adding the PDF and uncertainties in quadrature, as 



explained in Ref. [34 



The outline of the paper is as follows. Section [TTl reviews our NNLO implementation 
of DIS cross sections, placing an emphasis on the treatment of heavy-quark contributions. 
Section IIIII summarizes the statistical procedure of the global analysis, most notably the 
definition of the log-likelihood function and the implementation of correlated systematic 
errors. Section |IV] lists the experimental data used in the CTIONNLO analysis. Section 
IVl describes some features of the resulting CTIONNLO PDFs, while Section I VI I presents 
detailed comparisons of data and theory. Section IVlII presents NNLO predictions for collider 
measurements based on the CTIONNLO PDFs, and Section IVlIII contains our summary and 
conclusions. 



II. HEAVY-FLAVOR SCHEME IN THE CTIONNLO FIT 



A. QCD factorization for heavy quarks in the S-ACOT-x scheme 



A consistent implementation of contributions from the massive quarks (c and b) is a 
prerequisite and challenge for a viable NNLO PDF analysis. The mass dependence of the 
heavy-quark DIS contributions to PDF fits affects QCD precision observables in a wide 
range of energies [3|. Heavy-quark mass effects were studied in PDF fits since mid-1990's in 
the context of several theoretical approaches, or "heavy-quark schemes". The S-ACOTW; 
factorization scher ne l6| -l9| has been adopted in the recent NLO fits CTEQ6HQ [35|, 6.5 |3|, 



6.6 jl], and CTIO lo| . For the present work, the S-ACOT-x scheme has been extended to 
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NNLO accuracy, i.e., (9(a^), in the computation of neutral-current DIS cross sections 31 
The alternative TR' scheme [s^, 37 1 is used in the MSTW and HERAPDF fits, while the 
FONLL scheme has been adopted by the NNPDF collaboration. The BMSN scheme 



40, |4l| and the fixed-flavor number (FFN) schemes are used by the ABM and GJR groups, 



respectively. 

The general-mass variable flavor number (GM-VFN) schemes replace the zero-mass ap- 
proximation, which is no longer adequate for describing the DIS data. They evaluate the 
coefficient functions using the exact dependence on heavy-quark mass rrih, while heavy-quark 
PDFs provide an approximation for collinear production of hh pairs when ^ m^. They 
are valid across the whole range of Q values accessed in the global fits.^ 

Less variation in the predictions of various heavy-quark schemes was observed when going 
from NLO to NNLO calculations |3l|, H . For example, the spread of theoretical predictions 
for the standard- candle W and Z boson production cross sections at the LHC has shrunk 
from 6-8% at NLOj^ (with the largest deviations observed with the zero-mass PDFs) to 



2-3% at NNLO 




Quark- mass effects continue to be important. Their non-equivalent 



treatment by various PDF analysis groups gives rise to residual uncertainties in the standard 
candle predictions and in the observables sensitive to the gluon or heavy-quark scattering, 
such as Higgs boson production. 



The NNLO calculation in the S-ACOT-x scheme in Ref. (31| focused on two issues that 
had not been earlier addressed. First, we clarified the connection of the GM-VFN scheme at 
NNLO to the proof of QCD factorization for DIS with massive quarks presented by Collins 
[3]. It was demonstrated that the S-ACOT-x scheme can be derived to all orders from 
Collins' approach and is validated by a QCD factorization theorem. 

Second, we documented an algorithm that organizes the NNLO S-ACOT-x calculation in 
close analogy to the zero-mass VFN computation. As a result, the NNLO DIS cross sections 
can be constructed step-by-step from components that can be found in literature. In the 
S-ACOT-x scheme, all elements arise from the all-order factorization formalism, which is 
not always the case in other frameworks. These features distinguish the S-ACOT-x scheme 
from the TR' and FONLL schemes that follow different implementation approaches. 

As a result, the structure of the S-ACOT-x cross sections is readily reproducible; the 
universality of the PDFs follows from the Collins proof of QCD factorization The 
Q dependence of the QCD coupling and PDFs is found by numerical QCD evolution that 
assumes one shared A'^^ value in each Q range; QCD quantities for Nf and Nf + 1 active 
flavors are related at the switching energy scales through the matching conditions. The 
NNLO radiative contributions are assembled straightforwardly from the matrix elements in 
the massive FFN 17-1^ and zero-mass VFN schemes [H-ly], as well as from the mass- 
dependent operator matrix elements for the heavy-quark PDFs 



^ The chargcd-currcnt DIS cross sections, for which some 0(0^) massive amphtudes are not available, are 
evaluated at NLO. 

^ A complementary approach (an intermediate-mass scheme) uses approximate quark mass dependence in 
all scattering channels [42|. DIS cross sections in this approach utilzing the exact 0{as) massive ACOT 



terms and approximate 0{a^) and 0{a^) massive terms have been recently published 



43 
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The S-ACOT-x scheme reduces to the FFN scheme at the heavy-quark threshold ~ 
and to the zero- mass MS scheme at ^ m^, without additional renormalization. The 
matching of the GM-VFN cross sections to the FFN cross sections near the mass threshold 
is generally not automatic. To realize it, the TR' and FONLL scheme introduce additional 
elements (constraints on the Q dependence of the heavy-quark DIS contributions in the 
TR' scheme and the "damping factor" in the FONLL framework) that are not stipulated 
by the QCD factorization theorem. In the S-ACOT-x scheme, threshold matching is rather 
a consequence of the energy conservation condition that suppresses the difference between 
the GM-VFN and FFN results when approaches m^. This is achieved by restricting the 
allowed form for the approximate coefficient functions that describe QCD scattering off an 
initial-state heavy quark, so that they comply with energy conservation. It leads to the 
effective rescaling (x rescaling) of the light-cone momentum fraction in the approximate 
heavy-quark scattering terms 48| . 



The derivation of mass-dependent contributions from the QCD factorization theorem, 
which is shown to be compatible with the rescaled terms, leads to more confident constraints 
on the PDFs. Our fitting code implements these NNLO Wilson coefficient functions in the 
S-ACOT-x scheme together with the HOPPET program for the evolution of and PDFs 
m which the switching points between the active flavors can be expressed either in terms 
of the MS masses or pole masses.^ The fitting program can read either the pole masses or 
MS masses as an input. If the MS masses are selected as the input, they are converted into 
the pole masses when needed, e.g.^ to evaluate those operator matrix elements A^J^ that are 
published in terms of the pole masses. In the CTIO NNLO fit, the pole masses mc = 1.3 
GeV and = 4.75 GeV have been assumed. 



B. The S-ACOT-x scheme at NNLO in a nutshell 



In the S-ACOT-x scheme, a generic inclusive structure function F2^l = F takes form of a 
convolution product (®) of the Wilson coefficient functions Cia and the parton distribution 
functions fa/p{^,fi)- 



= Ee'E[c',,®/,/,](x,g), (1) 



where ^ is the light-cone momentum fraction, fi is the factorization scale, Nj indicates the 



number of active flavors, and Nf^ is the number of the produced final-state flavors (most 



The authors want to thank Gavin Salam for his help with the setup of this code. 
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generally, N^^ ^ Nj). The structure function can be also written as 



1=1 h=Ni+i 

where / and h are the indices of light-quark and heavy-quark flavors probed to the photon, 
respectively j39|. [Note that the S-ACOT-x scheme can simultaneously account for several 
flavors with non-zero masses]. On the right-hand side, 

F, = e2^[Q,,®/,/J(x,Q), F, = elJ2[Ch,a®fa/p]ix,Q). (3) 

a a 

The 0{al) radiative contributions, f/^'' and Fj:^\ are 

f/^) = ® if,, + fy,) + c-^'(^) ® E + eg ® , (4) 

{<f ® (A/P + k/p) + C^i? ® S + Ci^] ® , (5) 



F^'^ - el 



(2) 

where the lower-case notation c)^ ^ indicates a zero-mass Wilson coefficient function, 

the uppercase notation C^^l indicates a massive coefficient function, and S(a:,/i) = 

^fji [fi/p{x, fi) + fi/p{x, /i)] denotes the singlet-quark PDF. These equations have the same 
form as the factorized expressions for the zero-mass structure functions. Their components 



are listed explicitly in Ref. [31|. In this derivation, we employ a rescaling convention [9| to 

(k) 

construct the LO, NLO, and NNLO coefficient functions with initial-state heavy quarks, c\ ^ 
with k = 0,1, and 2, and the associated subtraction terms. They are obtained by evaluating 
the zero-mass expressions as a function of the rescaling variable x'- 

(k) fxQ mh\ (k) fx Q n\ at ^ c ^ i\ fa\ 

where 

X = x(^l + ^^j, (7) 

and Ylfs is the sum of the heavy-quark masses in the final state (equal to 2m/i in the 
lowest-order cc pair production). These rescaled coefficient functions obey energy conserva- 
tion and vanish near the production threshold, so that the FFN result is reproduced in this 
limit. 

Eqs. dl]) and ([5]) apply to the inclusive DIS functions F2, Fl, etc. In the case of semi- 
inclusive DIS production of heavy quarks, the definition of the semi-inclusive (SI) structure 
functions F/j 5/, such as F2'^ measured at HERA, requires additional care in order to obtain 
infrared-safe results at all O |4l| . In the CTIONNLO global fit, the following approximation 
for Fh^si has been adopted |3l|, l39|. At moderate Q values accessible at HERA, it is defined 
as 

Ni 

^S/(^, Q) = + E ^i^i''^'^ ® (/vp + fi/p)^ (8) 

1=1 
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where Ff^ {x,Q) is the 0{a1) contribution to the inclusive DIS function F{x,Q) arising 
from photon scattering on a heavy quark, cf. Eq. ([5]). The function L^^'^'^\^,Q / ^,mh/ ^i) is 
the non-singlet part of the light-quark component Fi{x,Q) that contains radiation of a hh 
pair in the final state, computed in Ref. [l9|. In the kinematic region Q < 10 GeV, which 

N S ('2') 

supports most of the HERA charm production data, we observe that Lj ^ contributes 
about 0-3% on the semi-inclusive charm cross section, i.e., it is small compared to the 
typical experimental errors. 



C. A phenomenological illustration 



In Fig. [H the CTIO NLO/NNLQ predictions for ^2^^ are compared to a recent data set 
from the HI collaboration S^], in which the charm structure function Fg'^ was extracted 



from the D* meson distribution. Predictions for F^'^ are computed as a function of the 
momentum fraction x for five different bins of Q. The dashed red lines and solid blue lines 
represent the S-ACOT-x predictions at NLO and NNLO, respectively. The inputs are the 
pole mass rric = 1.3 GeV and the strong coupling constant as{Mz) = 0.118. Although this 
specific data set (sol was not included in the CTIO fits, the overall agreement is very good. 
The NLO and NNLO predictions are close to one another for the most part, at least in the 
intermediate 10"^ < x < 10~^ kinematic region, in which the bulk of the data was collected. 
However, the NLO prediction is shown for a particular factorization scale that improved 
agreement with the data, n = ^/Q"^ + ml, while a more typical scale choice n = Q was 
taken at NNLO. In general, the NLO predictions for Fl'^ are characterized by a wide band of 
the scale dependence uncertainty. This band is reduced significantly when going to NNLO, 
and the NNLO predictions in the S-ACOT-x scheme are closer to those in the FONLL-C 
and TR' schemes [3l[. The NNLO predictions are therefore more robust compared to NLO. 



III. GLOBAL ANALYSIS WITH CORRELATED SYSTEMATIC ERRORS 



For completeness, and in order to establish notations for the ensuing discussion, in this 
section we summarize the statistical procedure adopted in the CTIONNLO analysis. An 
important aspect of this procedure is to determine the uncertainties of the PDFs, which 
arise from several sources, both experimental and theoretical. 

Experimental uncertainties may be uncorrelated between different measurements, such 
as the bin-by-bin statistical error on the measured cross section. Some experimental un- 
certainties are highly correlated between different measurements within one experiment or 
even among several experiments. The luminosity error is an example of a correlated exper- 
imental uncertainty in a data set. It affects equally the normalization of all cross section 
measurements from the experiment. Other correlated uncertainties typically exist and may 
have impact that is comparable to the luminosity error. 

Inclusion of the correlated experimental errors utilizing their standard deviations pub- 
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FIG. 1: The semi-inclusive heavy-quark function i^as a function of x in different bins of Q. The 
CTIONLO, CTIONNLO predictions and HI data Q are compared. 



lished by the experiments began in the CTEQ program with the construction of the CTEQ6 
parton distributions [H. This is achieved by constructing an appropriate figure-of- merit 
function that includes errors from both uncorrelated and correlated sources, as shown 
explicitly in Sec. 1111 A[ 

By examining the function in the neighborhood of the best fit and its dependence 
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on as, the CTIO global analysis determines a family of independent eigenvector sets that 
can be used to propagate the combined PDF+a^ uncertainty into theoretical predictions 
(cf. Sec. IIIIBp . The central PDF set and PDF eigenvector sets are generally sensitive 
to the implementation of correlated systematic errors in the fit. Possible procedures for 
reconstructing the correlation error matrix from the published standard deviations, and the 
expected modifications in the PDFs that they induce, are discussed in Sec. IIII CI A numerical 
comparison of these procedures will be presented in Sec. I VII 



A. Log-likelihood function in the CTIO NNLO analysis 



A typical experiment E in the global fit publishes a set of measurements {Mi;i = 
1,2,3, Npt}, consisting of a central value Di for the observable, a standard deviation 
for the uncorrelated experimental error Si, and standard deviations (3ka for each of Na sys- 
tematic errors, where k = 1,2,3, Na- We do not know the experimental errors, but 
experiments provide their standard deviations. So we write 



Di = Xi + Si6i + J2 PkaK (9) 



«=1 



where Xi is the "true" value of the observable.® Equation ([9]) defines nuisance parameters, 
6i and Aq,. These will be random numbers with mean value and standard deviation 1, 

(5,) = and = 1; (10) 
(A,) = and (A^) = 1. (11) 

We assume that the uncorrelated errors are completely uncorrelated among the measure- 
ments; that is, 

{6,6,) = 6i,. (12) 

We also assume that the systermatic errors are completely correlated among the Npt mea- 
surements; that is, Xa does not depend on i. However, the systematic arrors are among 
themselves uncorrelated, 

(AqA^) = 6ais. (13) 

The goal of the global analysis is to find the theoretical parameters for which the theoret- 
ical values Tj of the observables are as close as possible to the "true" values, X,. In Eq. 
replace X, by Tj. Then we are led to minimize the differences, by defining 

xl{{a},{X}) = xl + xl (14) 

where 

= 5^ ^ hOfc - Tfc -J2t^k,aK , (15) 

k=l k \ a=l / 



^ The "true" value is the mean value of Di that would result from a large number of independent experi- 
ments. 
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and 

and minimizing Xe with respect to the theory parameters. However, there are additional 
sources of systematic uncertainties associated with the unknown nuisance parameters {Aq}. 
Therefore, to have agreement between theory and data within the standard deviations of 
the experimental errors, we also vary the {Aq,} values, seeking to make Xe small. 

Because Xe only a quadratic function of {Aq,}, we may obtain the minimum of Xe with 
respect to {Aq,} analytically jsij. At the minimum, the nuisance parameters take the values 

Npt JVa 



where 



i=l 5=1 



A, = 5., + E^^. (18) 

k=l 

The corresponding best-fit representation for Xe is 

Npt 

minxl = ^(A - T,)(cov"^),,(Z}, - Tj). (19) 
It includes the inverse of the covariance matrix 



Nx 



(cov)ij = s'^5ij + A,a/3j,a, 



given by 



cov 



Nx 



(20) 



(21) 



The best-fit value of Xd ^^^^ expressed in terms of the data Dsh,i that are "shifted" 

from the central values by the best-fit correlated errors. 



N. 



pt 



i=l 



Here 



No, 



Dsh,i ^ ^ /^ia Aq, . 



(22) 



(23) 



(1=1 



we quote the function min x|) in Table [T] as a useful measure of the agreement between theory 
and data for each experiment. In Sec. IVit we will also use x\ = Xd + Xx another measure 
of the agreement between theory and data. 
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Thus far, we have considered only a single experiment. The global chi-square function 
that is minimized in the CTIO global analysis sums over all experiments E and includes a 
contribution Xthii^}) specifying theoretical conditions for the PDF parameters: 

xlioha.i = ^xl + Xth- (24) 

E 

The function Xth introduced to prevent some unconstrained PDF parameters from 
reaching values that might lead to unphysical predictions at small Bjorken x (specifi- 
cally, X < 10~^), where experimental constraints are sparse 0]. It rules out those PDF 
parameter combinations that may result in negative cross sections or unlikely flavor de- 
pendence. The specific condition imposed in the CTIO analysis is to constrain the ratio 
Rs{x, Q) = [s(x, Q) + s{x, Q)] I [^(x, Q) + d{x, Q)\ of the strange PDFs to non-strange sea 
PDFs to be in the interval 0.5 < Rs{x,Qq) < 1.5 at the initial scale Qq and x below 10~^. 
In the current fit, the Rg ratio is not constrained in this region by the experimental data, 
hence a loose theory-motivated constraint needs to be imposed. 

The global function in Eq. fl24|) is constructed by assuming that both uncorrelated 
and correlated errors are quasi-Gaussian and symmetric. When an experiment provides 
asymmetric errors, we symmetrize them. The exact procedure for symmetrization of errors 
has low impact on the outcome of the fit, as the number of points with very asymmetric 
errors typically is small compared to the total number of points and parameters. 



B. PDF eigenvector sets and ag uncertainty 



Besides the central (best-fit) PDF set, the CTIONNLO release includes PDF eigenvector 



sets to estimate the uncertainty range of our PDF fits using the Hessian method [5l|, |52 



and either the symmetric [1| or asymmetric [10|, |53| master formula to estimate the PDF 



uncertainties. The Hessian method is based on an iterative procedure for finding linear 
combinations of the fitting parameters by diagonalization of the Hessian matrix. In the 
CTIONNLO fit, we use 25 free parameters to describe the parton distributions at Qq, and 
hence have 25 eigenvector directions. For each of these directions, we find a pair of eigen- 
vector sets defined by moving away from the best-fit location (where takes its minimum 
value) by a distance that estimates the boundary of the 90 % confidence interval for each 
experiment. 



The method to find these eigenvectors is explained in [10|. In addition to including an 
upper "tolerance" bound on the increase in total that realizes the 90% c.l. agreement on 
average, as in the CTEQ6 analysis, we include a penalty term in x^ that quickly grows when 
the PDF set fails to describe any specific experiment. The effective function Xe// that is 
constructed this way is scanned along each eigenvector direction until xi// increases above 
the tolerance bound or quick xi// growth due to the penalty is triggered. 

The penalty term is constructed from statistical variables S derived from Xe values for 
individual experiments [lo[. In contrast to Xe? the variables S obey an approximate standard 
normal distribution independently of Npt, which simplifies the comparison of confidence 
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levels between data sets containing widely different numbers of points. The quasi-Gaussian 



S variables are found by using a simpler Fisher's approximation [5J] in the CTIO NLO 



fit (reliable for exp eriments with many data points, Npt ^ 9) and more accurate Lewis' 



approximation [55| in the CTIONNLO analysis. 



We note that even in the central fit, some data sets have xl;/Npt that lie outside the 90 % 
confidence level. That is not surprising, of course, since there are 28 data sets — naively the 
chance for all of them to lie with 90 % confidence is only 0.9^*^ = 0.05. To allow for this, for 
each data set that has Xe > ^pt the central fit, we rescale its Xe ^ factor Npt/x^entraifit 
before computing the penalty. 

In practical applications, the PDF uncertainty obtained with the Hessian eigenvector sets 
must be combined with the uncertainty. The procedure for combining these uncertainties 
is described in js^]. The central and error eigenvector sets of the CTIONNLO family assume 
as{Mz) = 0.118, which is compatible with the world-average value. In addition we provide 
best-fit PDF sets for other as{Mz) values in the interval 0.112-0.127. For a theory observable, 
the 90% c.l. as uncertainty can be estimated from the difference of predictions using PDF 
sets with as = 0.116 and 0.120. Then, the PDF-|-as uncertainty with all correlations can be 
estimated by adding the PDF and as uncertainties in quadrature. 



C. Implementation of correlated systematic errors 



In the treatment of correlated systematic errors, there is another subtlety concerned with 
the distinction between additive and multiplicative systematic uncertainties. The correlated 
systematic errors fall into two classes: additive errors, for which the experiment can deter- 
mine the absolute value Pi^a of the standard deviation; or multiplicative errors, for which 
only the relative fraction (Tj^Q = Pi^a/^i is known. The two kinds must be handled differently 
to avoid a bias in the outcome of the fit. However, it is a common practice in many exper- 
iments to publish the correlated systematic uncertainties as the relative percentage errors 
(Tj Q regardless of whether the systematic error is additive or multiplicative, rather than the 
absolute values f3i^a- 

To reconstruct the correlation matrix as 

Pi,a = 0'i,aXi, (25) 

one selects the reference central value Xj for each datum. Since the additive errors are 
supposedly independent of the theory predictions, it is natural, although not necessary, to 
use the experimental central values Di as the references, /3j_Q, = (Ji^aDi- 

For a multiplicative error, the reference to is generally unacceptable, as random fluctu- 
ations in Di tend to bias the best-fit parameters. A well-known multiplicative bias described 
by D'Agostini arises in the treatment of the normalization of the data {Hg], 13]. A fit that 



references the normalization error to Di would underestimate the true cross section. This 
downward bias is prevented if smoothly behaving Xi values are used as the references, such 
as the theoretical values Tj at each data point {Pi^a = (^i,aTi). 
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Most of the time the experimental paper does not distinguish between the additive and 
multiphcative errors. In this case, the global fit has to choose between several ways, and 
several tradeoffs, for computing /3j^Q,. When preparing the CTIO NNLO PDFs, we explored 
various procedures for the computation of f3i^a that have been identified in recent literature. 



[See also a related discussion in the appendix of Ref. [58 



1. Method D, or the "experimental normalization" method, normalizes all correlated 
errors to the experimental central values: that is, we compute (3i^a = o-i,aDi, for 
all additive and multiplicative errors alike. In this method, the correlated errors are 
independent of theory predictions. However, they are affected by irregular fluctuations 
of the central data points, which may result in a pronounced D'Agostini bias in the 
best-fit parameters when the fluctuations are large. 

2. Method T normalizes the multiplicative errors to the theoretical values that are up- 
dated in every fitting iteration. "The T'^^^ method" [s^ is a variation on this approach 
in which the multiplicative errors are updated once in many iterations. The T''^^ 
method was proposed to prevent nonlinear behavior of thai may occur in method 
T. Additive errors in these methods remain normalized to Di. Both methods T and 
T^^^ are free of D'Agostini's bias, as the multiplicative errors are estimated by using a 
smooth function. Method T*^"-* converges to the true solution after several updates of 
the t/''^ values, provided the partial derivatives dx^ /dTf^ are negligible, and may fail 
to do so otherwise. Our comparisons follow the original implementation of the T and 
y(o) methods j59|, in which only the luminosity errors (but not other multiplicative 
errors such as for the jet energy scale) were referenced to Tj {Tf^''). 

3. An extended version of method T normalizes both additive and multiplicative errors to 
the current theoretical values, /3j^Q, = ai,aTi. Similarly, the extended method T*^*^^ nor- 
malizes all correlated errors to fixed theoretical values, = ai^aT^^^ ■ The advantage 
of these methods is that Pi^a shows the smoothest behavior among all considered. 



The extended method T is used by default in our NLO fits. At NNLO the implementation 
of systematic effects became even more important, given the reduction in other uncertainties. 
In Sec. IVIDI we compare five procedures for implementation of systematic errors in inclusive 
jet production, where the systematic effects are among the most pronounced. We use the 
extended T method for other scattering processes. 



IV. EXPERIMENTAL DATA SETS AND THEORETICAL UPDATES 



The experimental data sets included in the CTIO NNLO fit are listed in Table HI With a 
small number of exceptions, they were chosen to be the same as in the CTIO NLO fit. At 
NLO, we presented two PDF sets, designated as CTIO and CTIOW. The distinction between 
them concerns the inclusion of the D0 Run-2 data for the rapidity asymmetry {Ai ) of the 
charged lepton from W boson decay. This data was included in the CTIOW analysis, with an 
extra weight, but not included in the CTIO analysis. More specifically, the CTIOW analysis 
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includes the Ag data points in three ranges of the lepton transverse momentum (pre) in the 
electron decay channel, and one pre bin in the muon decay channel. After the publication 
of the CTIO and CTIOW analyses, the D0 collaboration has recommended not to include 
those two less inclusive data sets in p^i, i.e., with 25 < pxi < 35 GeV and pti > 35 GeV. 
In the CTIONNLO analysis, only the most inclusive data with p^e > 25 GeV, in both the 
electron and muon decay channels, are included in the analysis. Since CTIONNLO includes 
only a part of the D0 Ai data that distinguishes between CTIO and CTIOW PDFs, it can 
be treated as an NNLO counterpart to either the CTIO or CTIOW NLO PDF sets. 

The Tevatron Run-2 inclusive jet data have a wider rapidity coverage (in the case of 
CDF) and smaller statistical and systematic errors, as compared to the Tevatron Run-1 
data. Their effects on the global analysis have been extensively discussed in Ref. [5|. We 
concluded that since the Run-1 jet data are not perfectly consistent with the Run-2 data sets, 
one may ask if the Run-1 jet contributions should be retained. In the CTIONNLO analysis, 
we do not include the Tevatron Run-1 inclusive jet data. This choice has implications for 
the large-x gluon PDF, as will be shown later. 

The data selection choices in regard to the Tevatron W asymmetry and inclusive jet 
production affect mainly certain combinations of the PDFs at a; > 0.1: specifically the 
ratio d/u (down- versus up-quark PDFs) and the gluon PDF. These combinations change 
by amounts that are mostly comparable to the PDF uncertainties obtained in the NLO 
analysis. 

Besides implementing the NNLO QCD contributions in neutral-current DIS and vector 
boson production, we updated the theoretical treatment of several experiments. In the 
coefficient functions for deep inelastic scattering, we have updated the definitions of the 
electroweak couplings describing Z and 7* — Z interference contributions. This update 
mildly modifies the d quark PDFs at large x. 

In the previous global analyses, the NLO jet cross sections were using tabulated point-by- 
point ratios of NLO/LO cross sections {K factors) computed with the EKS code |85|]. The 
original EKS calculation was published in the 1990's and had limited accuracy when compar- 
ing it to the latest precision jet data. A deeply modified version of the EKS code (MEKS) 
was prepared to provide advanced predictions for jet cross sections [s^. The MEKS cal- 



culation is entirely independent from NLOJET-I — |- [87|], another frequently used program 
for computation of NLO jet cross sections. A comparison to MEKS validated predictions 
from FastNLO and APPLGRID |9l|, the interfaces for fast interpolation of NLO- 



JET-I-+ cross sections. The specific input settings for which the MEKS, FastNLO, and 
APPLGRID programs agree to a few percent were documented 58, 8^. The magnitude of 



theoretical uncertainties in inclusive jet cross sections was also estimated. It was observed, 
for example, that the current NLO QCD scale dependence is too large at the highest Ptj 
and Uj of the LHC inclusive jet data set to provide meaningful constraints on the relevant 
PDFs. This situation is expected to be improved when the NNLO inclusive jet cross section 
calculation becomes available. 

In the CTIONNLO fit, we compute the NLO inclusive jet cross sections with FastNLO 
version 1.0 and cross- validate the results using MEKS. The two- loop threshold correction 
was not included. In the jet cross sections, both the factorization and renormalization 
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Experimental data set 


Npt 


CTIONNLO 


CTIOW 


Combined HERAl NC and CC DIS [60] 


579 


1.07 


1.17 


BCDMS F| [61] 


339 


1.16 


1.14 


BCDMS [62] 


251 


1.16 


1.12 


NMC F| [63] 


201 


1.66 


1.71 


NMC F^/F^ [63] 


123 


1.23 


1.28 


CDHSW F| [64] 


85 


0.83 


0.66 


CDHSW F^ [64] 


96 


0.81 


0.75 


CCFR F| [65] 


69 


0.98 


1.02 


CCFR xFaP [66] 


86 


0.40 


0.59 


NuTeV neutrino dimuon SIDIS [67] 


38 


0.78 


0.94 


NuTeV antineutrino dimuon SIDIS [67] 


33 


0.86 


0.91 


CCFR neutrino dimuon SIDIS [68] 


40 


1.20 


1.25 


CCFR antineutrino dimuon SIDIS [68] 


38 


0.70 


0.78 


HI F!^ [69] 


8 


1.17 


1.26 


HI 0-^ for CC [70, 71] 


10 


1.63 


1.54 


ZEUS F| [72] 


18 


0.74 


0.90 


ZEUS [73] 


27 


0.62 


0.76 


E605 Drell-Yan process, cr(pA) [74] 


119 


0.80 


0.81 


E866 Drell Yan process, a{pd)/{2a{pp)) [75] 


15 


0.65 


0.64 


E866 Drell-Yan process, cr{pp) [761 


184 


1.27 


1.21 


CDF Run-1 W charge asymmetry [77] 


11 


1.22 


1.24 


CDF Run-2 W charge asymmetry [78] 


11 


1.04 


1.02 


D0 Run-2 — )• efg charge asymmetry [791 


12 


2.17 


2.11 


D0 Run-2 VF. —7- charge asymmetry [80] 


9 


1.65 


1.49 


T "\ / X T r% ry ■ n-i i- i -i i- foil 

uv) Run-z Z rapidity distribution 


28 


0.56 


0.54 


CDF Run-2 Z rapidity distribution [82] 


29 


1.60 


1.44 


CDF Run-2 inclusive jet production [83] 


72 


1.42 


1.55 


D0 Run-2 inclusive jet production [841 


110 


1.04 


1.13 


Total: 


2641 


1.11 


1.13 



TABLE I: Experimental data sets examined in the CTIONNLO and CTIOW NLO analyses, to- 
gether with their values. 

scales are chosen to be equal to the transverse momentum Ptj of each individual jet. This 
choice is different from the one adopted in the CTEQ6.X and CTIO NLO analyses, where 
both scales were set to -PTj/2. 



V. THE CTIONNLO PARTON DISTRIBUTION FUNCTIONS 

With a similar setup as in the CTIO/CTIOW NLO analysis, the CTIONNLO fit results 
in about the same quality of the fit expressed in terms of the total chi-square . We obtain 
of order 2950 for 2641 data points, with minor variations dependent on the setup of the 
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fit. In the rest of tlie paper, we discuss representative results from the new fit, assuming 
NNLO PDFs, unless stated otherwise. 



X f(x,Q) versus x 





FIG. 2: CTIONNLO parton distribution functions. These figures show the Hessian error PDFs 
from the CTIONNLO analysis. Each graph shows 2: Uvaiencc = x{u — u),xdva.ience = x{d — d),0.10x g 
and 0.10 X Qsca Q-S functions of x for a fixed value of Q. The values of Q are 2, 3.16, 8, 85 GeV. 
The quark sea contribution is gsca = 2((i + u + s). The dashed curves are the central CTIO NLO 
fit. 



The central fit and PDF eigenvector sets are shown in Fig. [21 In addition to obtaining the 
central fit, there is also a need to understand the uncertainty ranges of the PDFs, resulting 
from the experimental uncertainties of the data included in the fit. For this purpose, we 
generate 50 alternate fits using the Hessian method. Figure [2] indicates the PDFs for the 
alternate fits, also called the "error PDFs", in addition to the central fit. In each graph, 
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PDFs for 4 parton combinations are shown, by plotting xf{x, Q) versus x for a fixed value 
of Q. The parton combinations are / = Wvaience, c^vaience, 9, Qsca.- The values of Q are 
2,3.16,8,85 GeV. The dashed curves are the corresponding functions for the CTIO NLO 
PDFs. Those can differ by significant amounts from the NNLO ones in some regions of the 
(x, Q) space, especially where x and Q are small. 




X 



FIG. 3: Ratios of various CTIONNLO central fit parton distributions to those of the CTIOW 
central fit, at Q = 2 GeV. 




FIG. 4: Same as Fig. H at Q = 85 GeV. 



To better assess these differences. Fig. |3] shows the ratios of various parton distributions 
from the CTIONNLO central fit to those from the CTIOW NLO central fit, at Q = 2 GeV. A 
few changes are noticeable: (1) At x < 10~^, the O(tts) evolution in CTIONNLO suppresses 
g{x,Q) and increases q{x,Q), as compared to CTIOW PDFs. (2) The heavy charm c{x,Q) 
and bottom b{x, Q) partons change as a result of adopting the 0{al) GM VFN scheme in 
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the CTIONNLO analysis. (3) In the large-x region, g{x,Q) and d{x,Q) are reduced by the 
removal of the Tevatron Run-1 inclusive jet data, the revised electroweak couplings, the 
alternative treatment of correlated systematic errors, and the scale choices. 

A similar comparison, but for Q = 85 GeV, is shown in Fig. |H It also indicates that at 
X > 0.1, d{x,Q), u{x,Q) and d{x,Q) are all reduced in CTIONNLO, compared to CTIOW. 
Likewise, g{x, Q) is also reduced in the large x region. 



g(x,Q) Q = 2GeV u(x,Q) Q = 2 GeV 




o 



o 



0.4 ' ' ' < < < < ' fcLKJ 

10"^ lO"* 10"^ 0.01 0.02 0.05 0.1 0.2 0.5 0.9 

X 



FIG. 5: Comparison of CTIONNLO error PDFs to CTIOW (NLO) error PDFs. All are expressed 
as ratios to the central fit of the CTIOW set. 



The uncertainties of CTIONNLO PDFs for various parton flavors are shown in Figure [3, 
compared to the CTIOW (NLO) uncertainties, at Q = 2 GeV. 



The CTIONNLO central fit is compared to that of MSTW2008NNLO (2^ in Figure Ei It 
shows that the CTIONNLO gluon and quarks are larger as x — )■ 0. At the smallest x values 
and at the initial scale Qq, the CTIO NNLO gluon remains positive, although consistent 
with zero, as a result of the chosen form of parametrization. We also note that CTIONNLO 
strangeness is larger at x ~ 10^'^. 
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FIG. 6: Ratios of various CTIONNLO central fit parton distributions to those of the 
MSTW2008NNLO central fit, at Q = 2 GeV. 
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VI. COMPARISONS TO INDIVIDUAL EXPERIMENTS 



As was already observed, the goodness of the NNLO fit to most experiments is about the 
same as for CTIO(W) NLO, and the changes in going to the NNLO analysis are comparable 
to the experimental uncertainties. The differences in Xe rarely exceed the expected statistical 
fluctuations of order ^y2Npt. Some improvement is observed in the fit to the HERA-1 
combined data on DIS, the Tevatron Run-2 jet production data, and the CCFR and NuTeV 
dimuon SIDIS. In this section, we present examples of the description of the data by the 
CTIO theoretical predictions, focusing on select precise measurements in DIS, W boson 
production, and jet production. 



HERA:e + P->e + X HERA : e + P -> e + X 

1.6, . . . . . . . 1 1.5| ^ ^ ^ ^ 




0.0 



1.0 1.5 2.0 3.0 5.0 7.0 1O0 1 50 20.0 2 5 10 20 50 100 200 

Q [GeV] Q [GeV] 

FIG. 7: The HERA combined data, with comparison to theory. The reduced cross section ar{x, Q) 
for NC DIS of positrons is plotted as a function of Q for 33 values of x. The graph on the left has 
X values in the range < x < 0.002; the one on the right has 0.002 < x < 0.65. The points are 
the published central values, while error bars correspond to the uncorrelated errors only. The red 
curves show the theoretical value of the reduced cross section aT.{x, Q) computed with CTIONNLO 
PDFs. 
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A. Deep-inelastic scattering at HERA 



The HERA combined data set for e^p neutral current (NC) and charged current (CC) 



deep-inelastic scattering (60| has 579 measurements of the reduced cross section ar{x,Q), 
after applying our restrictions Q > 2 GeV and W > 3.5 GeV, where the invariant mass of 
the hadronic final state is W = Q^/Xjx — 1. 

In Fig. [7] we show the e^p NC DIS reduced cross sections (366 data points) along with 
their uncorrelated bttot bars. Data are divided into two blocks: small x is shown on the left, 
i.e., < X < 0.002, and large x is shown on the right, i.e., 0.002 < a; < 0.65. Fig. [7] also 
shows theory curves based on the CTIONNLO parton distribution functions, superimposed 
on the the HERA combined data. Note that there are systematic differences between the 
central data and the theory, as there should be due to the systematic uncertainties of the 
data. Our detailed comparison of data and theory must indicate whether the systematic 
differences are consistent with the published errors [loj. 



HERA Data / Theory : e + P -> e + X 



4.0 r 
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0.0008 |-t-t» »«^>y«-t 



' 0.000574 



3 Q 0.0005 r^w^.M 



2.5 



2.0 



1.5 



1.0 



0.00041 



0.00032 ^ 



0.000246 



0.0002 
0.000161 



' 0.000124 



0.000105 



0.0000951 



0.000073 
0.000061 



10 



20 



Q [GeV] 



b 6 




5 10 20 50 100 200 
Q [GeV] 



FIG. 8: Ratio of data divided by theory, for the reduced cross section for positron-proton NC 
DIS. Theory is represented by the horizontal Unes. Indicated in both figures are the x-values: 
< X < 0.0016 in the left inset; 0.002 < x < 0.65 in the right inset. Black dots are experimental 
central values with uncorrelated errors only. Red dots are the optimally shifted data, where the 
shifts are evaluated according to the Hessian analysis of the systematic errors. 



Data divided by theory for e~^p NC DIS is shown in Figure [H in which we have a clearer 
illustration of the differences. The horizontal lines are marked by the corresponding values 
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of X. The black dots are the central experimental values, shown along with uncorrelated 
error bars. The red dots are the optimally shifted data values, determined by the Hessian 
analysis. We observe that the systematic shifts bring data and theory into better agreement, 
as expected. 

A histogram of the residuals of the shifted data, defined by 

Residuali = I A - ^ l^iaK " ) / = {Dsh,i - Ti) / Si, (26) 

is shown in Figure [91 The i-th residual is defined by the difference between the optimally 
shifted data and the theory, normalized to the uncorrelated error. The total number of 
systematic errors is Nx = 114. The curve shown in Fig. [9] is the standard normal distribution, 
with no adjustable parameters: it has a mean = 0, standard deviation = 1, and integral = 
1. The close similarity between the histogram of residuals and the Gaussian curve implies 
that the data and theory agree within the experimental errors, and that the experimental 
errors have been correctly estimated. 




FIG. 9: Histogram of residuals for the HERA combined data for positron-proton NC DIS. The 
residuals are defined in Eq. (j26p . The solid curve is the standard normal distribution, for compar- 
ison. 



The best-fit nuisance parameters {Aq,} for the HERA data are shown in the histogram of 
Figure [TOl The fact that the histogram agrees with the standard normal distribution shows 
that theory matches data within expectations based on the published correlated systematic 
errors. 
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Combined HERA DIS 



0.8 




Best-fit nuisance parameters 



FIG. 10: Histogram of the optimal normalized shift variables {Aq} for the HERA combined data. 
There are 114 systematic errors, all of which apply to all of the 579 measurements in the HERA 
combined data set. The solid curve the standard normal distribution, for comparison. 



B. W charge asymmetry and vector boson production 



The W charge asymmetry measured at the Tevatron is primarily sensitive to the shape of 
the ratio d{x, M\y)/u{x, M^) at x values above 0.1. The NNLO fit produces essentially the 
same quality of the agreement with the included Tevatron W asymmetry sets as the CTIOW 
NLO PDF set. The x^/^pt values for these data sets are listed in Table HI As an illustration, 
Fig. [m shows the CTIOW NLO and CTIO NNLO data residuals for the 12 points of the 



DO Run-2 lepton asymmetry data in the electron [79| and muon [80| decay channels. The 
residuals are very close for the two PDF sets. The quality of the NLO and NNLO fits to 
other data sets on vector boson production is also comparable. 



C. Inclusive jet production from Run-2 at the Tevatron 



The CTIONNLO global analysis only includes the Tevatron Run-2 inclusive jet data for 
the reasons explained in Sec. IIVI In what follows we examine the goodness of CTIONNLO 
fits to the Run-2 inclusive jet data. 

In Fig. [T2]we show the differential cross section d?a /{dyj dPTj) for inclusive jet production 



at D0 Tevatron Run-2 [8J|. The various curves correspond to the six rapidity intervals 
- 0.4, 0.4 - 0.8, 0.8 - 1.2, 1.2 - 1.6, 1.6 - 2.0, and 2.0 - 2.4. The points are the central data 
values, while the error bars (too small to be seen in most cases) represent uncorrelated errors. 
In the same figure, we also show the theory calculations based on the CTIONNLO parton 
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FIG. 11: Point-by-point residuals for the D0 Run-2 charge asymmetry data in the electron 
(upper inset) and muon (lower inset) channels. The residuals are defined in Eq. ()26p . 



distribution functions. In this figure we do not show systematic errors. With the systematic 
shifts, we obtain x|;/^pt = 115/110 for D0 Run-2, where 15 units are contributed by the 
penalty term x{ for the 23 nuisance parameters (Aq) controlling the systematic shifts, cf. 
the definition in Eq. (JH]). The overall effect of the correlated systematic shifts is significant, 
compared to the small uncorrelated errors, and is of the same order as in the CTIO NLO 
fits. More detailed comparisons, illustrating the impact of the systematic shifts on the jet 
data at NLO, can be found in jsj. 

To check the consistency of the comparison of the shifted data with CTIONNLO pre- 
dictions. Fig. [T3] shows histograms of the residuals, as defined in Eq. ( 126|) . for six bins of 
the jet rapidity {y). The distributions observed in the histograms appear to be reasonably 
consistent with the superimposed standard normal distribution, within the limitations of 
the small number of data points in each Uj bin. 

We also check the size of the optimal shift parameters {Aq;^ = 1, 2, 3, 23} and find 
that the histogram of the parameters shown in Fig. [14] is also consistent with the standard 
normal distribution. According to the figure, no unreasonable shifts are required to create 
agreement between the theory and the data. Some of the systematic errors have a larger 
impact on the jet cross section than others, such as those related to the jet energy scale 
uncertainties. 

A comparable measurement of inclusive jet production was performed by the CDF collab- 
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DO Run 2 inclusive jet data 

100 , 
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PT[GeV] 

FIG. 12: Data and theory for the inclusive jet production data at D0 Tevatron Run-2. The black 
points represent central values of the data, error bars are uncorrelated errors, and red curves are 
NLO theory obtained by using CTIONNLO PDFs. 

oration in Run-2 j83|. This jet measurement has different systematic uncertainties than D0, 
and hence provides an independent constraint on the PDF parameters. A similar investiga- 
tion indicates that the CTIO PDFs agree well with the CDF measurement. Fig. [T5] compares 
the CDF data set for inclusive jet production and the CTIO prediction, and the correspond- 
ing histogram of the residuals is shown in Fig. [161 The CDF value of x\/^pt = 102/72 = 1.42 
is higher than for D0, but the histograms of the residuals do not reveal pronounced system- 
atic disagreements. The contribution of 25 systematic shifts to the xl; function of the CDF 
measurement is Xa = 18- 



D. Correlated systematic errors in Tevatron jet production 

The results of the previous subsection show the importance of experimental systematic 
uncertainties in describing the Tevatron jet data. These uncertainties can be generally 
of two kinds, additive and multiplicative. For example, the underlying event and pileup 
uncertainties are additive errors. They are determined independently of the jet cross section. 
The luminosity uncertainty and jet energy scale (JES) uncertainty are multiplicative errors. 
From the start, they are estimated as the percentage of the jet cross section in each bin. In 
jet production, there is no reason to handle the luminosity uncertainty differently from the 
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FIG. 13: Histograms of residuals for the D0 Tevatron Run-2 inclusive jet production data. Resid- 
uals are defined in Eq. (|26p . The solid curve is the standard normal distribution, for comparison. 



JES and other multiplicative errors, as most of them have a comparable impact on the fit.^ 
Often, it is unknown if the error is additive or multiplicative. All of them are reconstructed 
from the published percentage values according to one of the methods that were reviewed 
in Sec. ITTTn 

Both the QCD scale dependence of NLO jet cross sections and the detailed implemen- 
tation of the systematic errors were explored when preparing our NNLO set of PDFs. An 
assortment of NNLO gluon PDF shapes, obtained by various treatments of jet production 
data in the global analysis and normalized to the shape of the CTIOW NLO gluon PDF, is 
illustrated in Fig. [T71 We compare five methods for the computation of Pi^a in jet production 
(cf. Sec. IIII Cp and several choices of the renormalization and factorization scale. The CTIO 



^ In the CTEQ6.X series and some other previous fits, the normalization of each experiment was fitted 
as a separate free parameter with a penalty imposed on its deviation from the nominal value. The 
procedure used was similar to the "normalization penalty trick" discussed in Ref. [s^. In the CTIO fits, 
the normalization uncertainty is instead included in the correlation matrix /3i_a and treated similarly to 
other multiplicative errors. 
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DO Run-2 inclusive jet production 




Best-fit nuisance parameters 



FIG. 14: Histogram of the optimal systematic shifts {Aq} for the D0 Tevatron Ruii-2 inclusive 
jet production data. 

central set and its PDF error band, both obtained by assuming method D for computing 
fii^a and the QCD scale /i = P^j, are shown by the solid black line and shaded green region, 
respectively. For data sets other than the jet data in the global fit, we consistently use the 
extended T definition of 

If the Run-2 jet cross sections are not included in the calculation of the global (o^' 
included with a small weight, such as 0.01), the fit generally prefers a softer gluon PDF at 
X > 0.1, as indicated by the magenta long-dashed line in the figure. If as{Mz) is fitted as 
well, its NNLO value in a fit without the jet data sets is also low, of order 0.114. 

With the jet data included, there is a general preference for an enhanced g{x, Q) at large x, 
with some possible outcomes shown in Fig. [171 The curves were obtained using the extended 
T, T, and D methods, for the factorization scales /i = Ptj/^, PTj,2PTj, with Ptj being the 
transverse momentum of an individual jet in each Pt bin. Both the original method T 
(applied only to the luminosity error) and extended method T (applied to all correlated 
errors) produce a robust increase in g{x, Q) at x > 0.2, but the shape and magnitude of the 
enhancement vary with fi. The D method (black solid curve) results in a softer shape of the 
best-fit g{x,Q) which lies between the fits with and without the jet data. All these curves 
lie well within the PDF error band of the public CTIO NNLO set. 

A complementary perspective is provided by Fig. [181 comparing the gluon PDF g{x, Q) 
at Q = 85 GeV obtained with the five methods for computation of (3ia and the same 
factorization scale /i = Ptj in FastNLO jet cross sections. All the curves are normalized 
to the CTIO NNLO prediction (based on the D method). Here, we observe some differences 
between the (extended) T and T^^^ prescriptions. The extended T method results in a 
drastic enhancement of g{x, Q) at x > 0.3, which is tempered when using the extended T^^^ 
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FIG. 15: Data and theory for inclusive jet production at CDF Tevatron Run-2. The black points are 
the central values of the data, the error bars are uncorrelated errors, and the red curves represent 
the NLO predictions computed by using CTIONNLO PDFs. 



method. The T and T^^^ methods produce a stronger enhancement at 0.05 < x < 0.3, but 
are flatter for x > 0.3. 



Parameter 


Type of 
the fit 


Method for computing f3i^a 
D T r(0) Ext. T Ext. r(°) 


Xl/Npt for CDF {Npt = 72) 
and D0 {Npt = 110) 


PDF only 
PDF+Q, 


1.42, 1.04 1.23, 1.05 1.21, 1.04 1.48, 1.03 1.49, 1.03 
1.33, 0.95 1.24, 1.07 1.23, 1.06 1.48, 1.06 1.51,1.08 


Alumi 


PDF only 
PDF+Q, 


1.2,1.3 -1.4, -1.2 -1.5, -1.3 0.53, 0.90 0.61, 0.96 
2.7,2.8 -0.71, -0.51 -0.75, -0.55 1.2,1.5 1.4,1.8 




PDF only 
PDF+Q, 


0.118 0.118 0.118 0.118 0.118 
0.1137 0.1159 0.1159 0.1162 0.1159 



TABLE II: Best-fit values of XE/^pt ™d luminosity nuisance parameter Aiumi for the Tevatron 
Run-2 inclusive jet production, and for as{mz) in a fit with floating a<j(M^). The columns are 
obtained with the five methods for including correlated systematic errors in jet production. The 
first (second) numbers in each column correspond to the CDF (D0) data sets, respectively. 



In all of the listed fits to the jet data, the quality of the fits are generally comparable. 
Table [ITl shows representative output quantities for five methods: /Npt for the CDF Run-2 
and D0 Run-2 jet data sets, the best-fit luminosity parameter Aiumi? and the best-fit as{Mz) 
(if fitted). In all fits, the values of /Npt are similar, with the T/T^'^^ methods showing 
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FIG. 16: Histogram of residuals (defined in Eq. (j26p ) for the CDF Run-2 inclusive jet production. 



a marginally better agreement with the CDF jet data {x^/Npt ~ 1.23 vs. 1.4-1.5). If the 
QCD coupling strength is fixed in the fit to be as{Mz) = 0.118, all luminosity shifts do not 
exceed 1.5 standard deviations, see the lines referring to "PDF only" fits. 

However, if as{Mz) is also fitted (the 'TDF+a^" fits), in the D method we obtain a low 
as{Mz) of 0.1137, accompanied by a large luminosity shift Aiumi ~ 2.7. The D method also 
results in a low x^/^pt combination of 1.33 and 0.95. 

With four other methods, we obtain higher best-fit values of as{Mz) ~ 0.116, higher 
values, and the luminosity shifts below 2a. This peculiarity of the "FDF+a^" fit employing 
the D method is suggestive of a strong D'Agostini bias. It does not occur at the same level 
in the "PDF only" fits that use a fixed a,(M^), including CTIO NNLO. 

As can be seen from Fig. [171 by using method D to determine the published CTIO NNLO 
set and the QCD scale /i = Ptj, we obtain a PDF uncertainty band that also covers the 
gluon PDF parametrizations found with the alternative methods for the correlation matrix 
and choices of QCD scales. The softer gluon behavior than in CTIOW NLO is consistent 
with the preference of the non-jet experiments and results in lower for the combined 
HERA DIS set. It may also simulate the anticipated effect of missing NNLO corrections 
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FIG. 17: Comparison of NNLO gluon PDF parametrizations obtained using different procedures 
for including jet production data sets in the global fit. 

g(x, Q = 85 GeV) normalized to CT10 NNLO 
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FIG. 18: Dependence of the gluon PDF on various definitions of the correlation matrices Pia for 
the Tevatron Run 2 jet data in the global fit. 

on the jet cross sections. The uncertainty band for g{x, Q) in Fig. [17] thus represents the 
current level of uncertainty in predicting the inclusive jet cross sections, due to variations 
in the fit procedure. These variations are typically smaller than the PDF error determined 
by the Hessian method. 



31 



VII. COLLIDER PHENOMENOLOGICAL PREDICTIONS 



A number of standard model cross section measurements from the LHC have been pub- 
hshed, at center-of-mass energies of 7 and 8 TeV. Differential distributions are available 
from the 7 TeV data (taken in 2010 and 2011), while the current 8 TeV results (2012) are 
mostly in the form of inclusive cross section measurements. In this section, we use the CTIO 
PDFs presented in this paper to determine theory predictions for several important LHC 
cross sections, including the production of electroweak gauge bosons, inclusive jets, and top 
pairs. We compare the resulting predictions with the most recent corresponding LHC data 
sets. These data sets have not been used in the determination of the CTIO PDFs, but will 
be used in future fits. 



We also compare some of the CTIO predictions to those based on other NNLO PDFs: 

All calculations use the native 



NNPDF2.3 [25[ 
value of as(mz 



HERAPDF1.5 28,29 



and MSTW2008 [2J 
for the PDF set; the PDF uncertainties are shown at the 68% c.l. For the 
results from HERA PDFs we only include the uncertainties from the standard eigenvector 
sets. We do not consider the a, uncertainties. 



The LHC standard model cross sections are calculated with the following settings. 



For electroweak gauge boson production, 
NNLO using the program FEWZ2.1 [oi. 
The differential cross sections and asymmetries with cuts on transverse momenta of 



the total cross sections are computed at 
94l | and the factorization scale Q = My- 



96 that realizes the NNLL re- 



the leptons are computed with the code ResBos [95 . 
summation of logarithms at small transverse momentum. The settings of the ResBos 
computation and nonperturbatitve function are the same as in [93] , including the hard 
QCD scale that is equal to Mv/2. This calculation provides a close approximation 
to the NNLO resummed cross section. The needed two-loop coefficients are obtained 



from the 0{a'^) cross sections at large transverse momentum |98|, |99| a nd b y requiring 



that the ResBos tot al rate co incides with the inclusive NNLO rate [100[ | computed 



by the Candia code [10l|, llQ2 



The cross se ction for top quark pair production is computed by the TOP-I-+ code 



version 1.4) |l03 |. The factorization and renormalization scale is chosen as Q = rrit 



the top quark mass is = 172.5 GeV |l04 |. and the other settings are as in Ref. |l05 



The calculation includes the exact NNLO corrections to the qq — >■ it production |106l 
and an approximate NNLO combined with NNLL threshold logarithms in the other 
channels. Hence, the calculation provides a close to exact NNLO prediction at the 
Tevatron, where the qq channel dominates, and an approximate NNLO prediction at 
the LHC, where top production takes place primarily from gg initial states. 

Higgs boson production cross sectio ns in the g luon-gluon fusion channel are computed 



at NNLO using the HNNLO code |107l . Il08| . The factorization and renormalization 
scale was taken to be Q = fnu = 125 Ge V, co nsistent with the recommendation of the 



LHC Higgs cross section working group fiogl . Higgs production in the bottom quark 
fusion channel is computed at NNLO using the bbh@nnlo vl.3 code liol . 
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LHC inclusive jet cross sections are computed at NLO with the FastNLO code 88 

and cross checked against the MEKS code in [58|, |86|. The central scale in these 
predictions is chosen to be /i = Ptj, the transverse momentum of each jet. 
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FIG. 19: NNLO W+ and W cross sections at the LHC. 



In Fig. [ini we show prediction s at NNLO for production at 7, 8 and 14 TeV, com- 
pared to ATLAS data at 7 TeV [ml and CMS data at 8 TeV [ll2|. The CTIO NNLO 
predictions are in good agreement both with the data and predictions from the other NNLO 
PDFs listed. The central values of the data are below all NNLO predictions, but are within 
the uncertainties. Among all PDF sets, HERAPDF predicts the largest cross sections. For 
simplicity, the error bars for the theoretical predictions indicate the symmetric PDF uncer- 
tainties at 68% c.l. Generally, the PDF errors of the and W~ are strongly correlated 
0], but the 68% c.l. error ellipses overlap well for all PDF sets. 



In Figure [201 we show the NNLO predictions for the W and Z total cross sections at 
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FIG. 20: NNLO and Z cross sections at the LHC. 



7, 8 and 14 TeV, again compared to ATLAS data at 7 TeV and CMS data at 8 TeV 



112| . The conclusions are broadly similar to those for the and W comparisons. 



B. Top quark pair production, total cross sections 



In Fig. fITi we sho w the com bined Tevatron total cross section for tt production |113 



and ATLAS |114| - |118| and CMS |119j total cross sections, compared to approximate NNLO 
predictions. The experimental central values and uncertainty ranges are indicated by vertical 
lines and colored rectangles. For each collider energy, all cross sections are normalized to the 
respective CTIONNLO prediction. Again, we observe good agreement between CTIONNLO 
and other theoretical predictions as well as the data. 
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FIG. 21: Approximate NNLO predictions and experimental measurements of tt total cross sections 
at the Tevatron and LHC, normalized to CTIONNLO predictions. 



C. W and Z rapidity distributions 



More detailed PDF information can be gained by c ompa risons to W /Z rapidity distri- 
butions and W charge lepton asymmetry from ATLAS [m|, CMS [l20|, and LHC-B [m 
at 7 TeV. For differential observables, we use ResBos to account for multiple soft gluon 
radiation that may be non-negligible in the presence of constraints on the charged lepton's 
transverse momentum pxe- Typical W/Z production measurements require ptc to be above 
20-25 GeV in order to suppress charged leptons from background processes, and because of 
trigger requirements. 

In some measurements, the lower pt£ cut is as high as 35 GeV, approaching the position 
of the Jacobian peak in da/dpn at pxi ~ M]y/2 ^ 40 GeV (for W boson production), 
where large logarithms ln(Qr/Q) dominate. {Q and Qt denote the invariant mass and 
transverse momentum of the weak boson, respectively.) With such high pxe cuts, the Res- 
Bos resummation calculation, which resums large logarithmic contributions to all orders in 
strong coupling, is expect ed to provide a better description of the data than the fixed-order 



calculations |93l . Il22l . Il23l |. In practical comparisons, differences between the resummed and 
fixed-order NNLO calculations are most pronounced with an upper cut close to pti ~ (5/2 
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FIG. 22: Rapidity distributions of Z and cross sections in ATLAS (35 pb ^) measurements 
at 7 TeV. 



124j . The resummed and NNLO asymmetries are closer in general when only a lower cut 



is imposed. 



In Fig. [22] we show a comparison between predictions for the Z and lepton rapidity 
distributions and the data measured by the ATLAS collaboration using 35 pb~^ of integrated 



luminosity at 7 TeV . In all three figures the green band represents the PDF uncertainty, 
while the red solid line shows the central prediction. For the Z rapidity distribution, the 
central ResBos prediction overshoots the data by a few percent, but is within the PDF 
uncertainties from the data. These LHC Z/W± data are yet to be included in the fits and 
may modify the quark PDFs at small x. It is interesting to note that ATLAS found that the 
tight constraints on the relevant normalization between the W and Z data sets resulted in 
a larger-than-expected strang e qu ark density in the x range represented by the Z and 



lepton rapidity distributions |125| . We shall include these (and upcoming) data sets, with 



correlated systematic error analysis, in our next run of global analysis. 



In Fig. [231 a Res Bos prediction is compared to the lepton charge asymmetry from 
ATLAS and LHC-B |l21| at 7 TeV as a function of the lepton pseudorapidity r]iep. The 
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FIG. 24: A lepton rapidity distribution of the W charge asymmetry in the CMS (840 pb ) 
measurement at 7 TeV. 
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FIG. 25: PDF-induced correlations between W charge asymmetry Achiiliep) at the LHC 7 TeV and 
Uy{x,Q) = u{x,Q) — u{x,Q), dy{x,Q) = d{x,Q) — d{x,Q, and d{x,Q)/u{x,Q) at Q = 85 GeV. 
The charge asymmetry is evaluated at lepton rapidity values indicated in the figure. 

agreement of CTIO NNLO with both asymmetry sets is good, e.g. Xs/^pt = 0.45 for 
ATLAS W asymmetry. The LHC-B experiment probes a different kinematic range (more 
forward coverage) than either ATLAS or CMS, and thus can provide unique information 
on large-x quark distributions. In this figure, the central prediction and PDF error are 
indicated by the black solid and blue solid lines, respectively. 

The LHC W charge asymmetry measurements have the potential to test some combi- 
nations of the PDFs that are still poorly constrained. For example, the measurement of 
the charge asymmetry in the rapidity bins of CMS in Fig. probes valence quark PDFs, 
Uy{x, Q) = u{x, Q) — u{x, Q) and Q) = d{x, Q) — d{x, Q) in the small-x region that is 
not constrained by the previous data. 

Fig. 123 shows the PDF-induced correlation cosine cos ^p 0, Ell between the LHC W 
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charge asymmetry in several bins of lepton rapidity rjiep and valence u and d quark PDFs 
and the ratio d{x, Q)/u{x, Q) of the CTIONNLO PDF set. In the intervals of x where cosyj 
is close to ±1, the W asymmetry data at the shown rjiep values is sensitive to the PDF 
uncertainty on the indicated parton distributions. A strong correlation can be observed 
with Ui, at 0.01 < x < 0.15 and at 0.005 < x < 0.05. Both valence PDFs show a strong 
anticorrelation at x < 0.001 as a consequence of the valence sum rules. The d/u ratio shows 
a strong correlation at 0.0005 < x < 0.005 and anticorrelation at 0.05 < x < 0.2. The other 
PDF flavors are not strongly correlated with the W asymmetry in these rapidity ranges. 

The current CMS W asymmetry data shown in the left panel of Fig. agrees reason- 
ably with the CTIO NNLO PDFs, given the present PDF uncertainties and experimental 
uncertainties. The CMS data is slightly below the central CTIO prediction at \r]iep\ < 1, but 
within the PDF uncertainty band. The disagreement is increased at |?7/ep| < 2, where we 
also observe less regular behavior of the data. 

We note that the reconstruction of the present CMS data has a systematic uncertainty 
associated with the discrimination from background processes. In the Tevatron W asym- 
metry measurements, it was essential to require a sufficiently high missing (MET), of 
order 20 GeV or more, to discriminate the W boson decay events from the significant back- 
ground. The MET requirement has not been imposed in the current CMS measurement, 
which applied a different technique to suppress the background. 

In the ResBos theoretical calculation, we cannot exactly implement the experimental 
background subtraction technique when comparing to the CMS data. At the lowest order in 
PQCD (when Qt = 0), the additional MET cut is not needed, as the condition pxi > 35 GeV 
automatically implies MET > 35 GeV. At higher orders, the pxi and MET cuts are no longer 
equivalent. The predictions for W charge asymmetry depend on the assumed MET cut, cf. 
the right panel in Fig. |2H This behavior suggests that subtle effects in the separation of W 
boson events from the background may be comparable to the observed differences between 
CTIO NNLO theory and CMS data in Fig. . 

In the absence of applicable experimental constraints in the CTIO fit, the small- a; behavior 
of the and d^ PDFs is assumed to be governed by a shared x power, m„(x, Qq), d^^x, Qo) ~ 
^YiQ X with the same Ai^ ^ This Regge-inspired assumption partly 

explains the correlations of PDFs at small x. It will be tested by including the LHC W 
asymmetry in the upcoming fits. 



D. ATLAS inclusive jet distribution 

In Figs. [26] and |271 we compare the inclusive cross section predictions to the ATLAS 
measurements using the anti-Zcr jet algorithm with R = 0.6. In Fig. |26l the comparisons 
are made to the ATLAS raw, unshifted data; in Figs. [271 the optimal systematic shifts have 
been applied to the ATLAS central values. The scale uncertainty in the theory prediction 
is calculated by varying the QCD scale in the range -Prj/2 < < "^Ptj- The uncertainties 
of the data points shown in Figures [21] and [27] are evaluated by adding the statistical and 
uncorrected systematic errors in quadrature. 



39 
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FIG. 26: Comparison of CTIONNLO predictions with unshifted 2010 ATLAS inclusive jet data 
{R = 0.6). 




200 400 600 800 1000 1200 100 200 300 400 500 



Pt (GeV) Pr (GeV) 



FIG. 27: Comparison of CTIONNLO predictions with shifted 2010 ATLAS inclusive jet data 
{R = 0.6). 



The overall agreement of the shifted data with the (NLO) predictions utilizing 
CTIONNLO PDFs is good even without theory uncertainty: xl/Npt = 0.74 (0.78) for 
R = 0.6 (0.4). Some differences in the shape and normalization observed for the unshifted 
data disappear upon application of the systematic error shifts. The Xe function is defined 
according to the extended T convention for the ATLAS correlation matrix. Slightly higher 
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values of Xe would be obtained with alternative conventions for /3j Q, 58|. An in-depth com- 
parison of theoretical jet cross sections based on CTIONNLO to the 2010 ATLAS data and 
predictions based on other PDF sets is documented in jHsf. In that paper, we examine the 
uncertainties in the current NLO predictions compared to the experimental accuracy and 
discuss prospects for constraining the PDFs with the upcoming LHC jet data. 



E. SM Higgs boson total cross sections 
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FIG. 28: NNLO cross sections of SM Higgs boson production via gluon fusion (left) and bottom 
quark annihilation (right), normalized to CTIONNLO predictions. 



Now we turn to the NNLO total cross se ctions for Higgs boson product ion a t both the 



107 , 108 1 and bb annihilation |llO |. shown as 



Tevatron and LHC via gluon-gluon fusion 
ratios to the CTIONNLO prediction in Fig. |28l At the Tevatron, PDF uncertainties of 
the Higgs boson production cross section via gluon-gluon fusion are large and cover all the 
central predictions. At the LHC, the CTIONNLO gluon-gluon fusion predictions are in a 
good agreement with those from the other PDF groups, except for NNPDF2.3 at 7 and 
8 TeV, which predicts larger cross sections. For Higgs boson production via bb annihila- 
tion, the NNLO predictions obtained with all PDFs are consistent within the quoted PDF 
uncertainties. 



VIII. DISCUSSION AND CONCLUSION 



We have presented next-to-next-to-leading order (NNLO) parton distribution functions 
(PDFs) from the CTEQ-TEA group. These CTIONNLO PDFs have been determined based 
on essentially the same global data sets used in the previous CTIO and CTIOW NLO PDF 
analyses. In this new analysis, the effects of finite quark masses have been implemented in 
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the S-ACOT-x scheme at NNLO accuracy. We obtain a similar quahty of agreement with 
the fitted experimental data sets in the NNLO fit as at NLO. 

We find that at low x (below 10~^), the NNLO gluon distribution is suppressed, while the 
quark distributions increase, compared to the same distributions at NLO. The 0{al) GM 
VFN scheme used in the NNLO fit results in changes to the heavy quark distributions, both 
charm and bottom. The large- a; gluon and ci-quark distributions are reduced due to (1) the 
removal of the Tevatron Run-1 inclusive jet data, (2) the alternate treatment of correlated 
systematic errors and choices of renormalization and factorization scales in jet cross sections, 
and (3) revised electroweak couplings in DIS cross sections. 

Compared to the MSTW2008 NNLO PDFs, the gluon and quark distributions are larger 
as X approaches zero (the CTIO parameterization requires a positive gluon), while the 
strangeness distribution is larger over most of the x range. The differences tend to decrease 
as increases. 

We have compared NNLO predictions using CTIO to available LHC data and have found 
good agreement. Numerical comparisons presented here complement a more detailed bench- 
mark study of the LHC predictions that has been recently released j58|. The available LHC 
data is already providing important information on PDFs, and future data will provide even 
stronger constraints. This will be developed further in a future publication. 
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